CoCalc Logo Icon
DocsShareSupport News Sign UpSign In
Views: 100
Image: ubuntu2004
Embed | Download | Raw |
Kernel: Python 3 (ipykernel)
import random
from IPython.core.display import SVG
import pyomo.environ as p
from pysat.solvers import Solver
from pysat.formula import CNF
import py_svg_combinatorics as psc
from ipywidgets import widgets, HBox
from collections import Counter
from pprint import pprint
from random import randint
import numpy as np
from IPython.display import IFrame
import IPython
from copy import copy
import os
from pathlib import Path
from collections import defaultdict
from itertools import combinations
import time


from IPython.display import IFrame
title_ = 'Minimum Test Collection'
IFrame(f'https://discopal.ispras.ru/index.php?title=Hardprob/{title_}&useskin=cleanmonobook', width=1280, height=700)

model = p.ConcreteModel()

Визуализация

Без ограничения общности, пусть у нас подмножество SS – численное.

model.InitialSet = range(10)
model.SetLen = len(model.InitialSet)
model.SetIndex = range(model.SetLen)

Зададим семейство подмножеств СС:

model.InitialSubsets = [[1,2], [1,2,3], [5,6], [6,8,], [0,4,7], [0,1,2,3,4,5,6,7,8,9], [0], [1], [2], [5,7]]

svg = psc.subsets2svg(model.InitialSubsets)
SVG(data=svg)
Image in a Jupyter notebook

Запомним его мощность и проиндексируем подмножества.

model.SubsetCount = len(model.InitialSubsets)
model.SubsetIndex = range(model.SubsetCount)

Удобнее будет хранить его в виде таблицы связности:

S[i][j]=1S[i][j] = 1, если ii-ое подмножество содержит элемент jj

S = [[0 for _ in range(model.SetLen)] for _ in range(model.SubsetCount)]

for i in range(model.SubsetCount):
    for elem in model.InitialSubsets[i]:
        S[i][elem] = 1

model.Subsets = S

Реализация в Pyomo

Переменная x[i]=1x[i] = 1, если мы берем ii подмножество в целевое семейство C′C' и 00 в противном случае.

model.x = p.Var(model.SubsetIndex, within=p.Binary)

def display_model(model):
    selected_sets = [subset for subset in model.SubsetIndex if round(model.x[subset]()) == 1]
    return SVG(psc.subsets2svg(model.InitialSubsets, selected=selected_sets))

Целевая функция — мощность семейства C′C', т.е. сумма x[i]x[i] по ii

@model.Objective()
def FinalCapacity(model):
    return sum(model.x[subset] for subset in model.SubsetIndex)

Единственное ограничение:

∀i,j∈S  ∃Ck∈C′:[i∈Ck,j∉Ck,i∉Ck,j∈Ck\forall i,j \in S\; \exists C_k \in C': \left[ \begin{gathered} i \in C_k, j \notin C_k, \\ i \notin C_k, j \in C_k \\ \end{gathered} \right.

То есть

∀i,j∈S∃Ck∈C′:  Subsets[Ck][i]XOR  Subsets[Ck][j]==1\forall i,j \in S \exists C_k \in C':\; Subsets\left[C_k\right]\left[i\right] XOR \;Subsets\left[C_k\right]\left[j\right] == 1

Теперь можно просто заменить существование на сравнение суммы по k с единицей.

@model.Constraint(model.SetIndex, model.SetIndex)
def all_pair_of_elements_are_included_excluded(model, i, j):
    if i <= j:
        return p.Constraint.Feasible
    else:
        return sum(model.x[subset] * (model.Subsets[subset][i] ^ model.Subsets[subset][j]) for subset in model.SubsetIndex) >= 1

Определяем солвер:

solver = p.SolverFactory('cbc')
solver.solve(model).write()
# ========================================================== # = Solver Results = # ========================================================== # ---------------------------------------------------------- # Problem Information # ---------------------------------------------------------- Problem: - Name: unknown Lower bound: 8.0 Upper bound: 8.0 Number of objectives: 1 Number of constraints: 3 Number of variables: 3 Number of binary variables: 10 Number of integer variables: 10 Number of nonzeros: 3 Sense: minimize # ---------------------------------------------------------- # Solver Information # ---------------------------------------------------------- Solver: - Status: ok User time: -1.0 System time: 0.01 Wallclock time: 0.01 Termination condition: optimal Termination message: Model was solved to optimality (subject to tolerances), and an optimal solution is available. Statistics: Branch and bound: Number of bounded subproblems: 0 Number of created subproblems: 0 Black box: Number of iterations: 0 Error rc: 0 Time: 0.49807310104370117 # ---------------------------------------------------------- # Solution Information # ---------------------------------------------------------- Solution: - number of solutions: 0 number of solutions displayed: 0 
selected_sets = [subset for subset in model.SubsetIndex if round(model.x[subset]()) == 1]
selected_sets
[0, 1, 2, 3, 4, 6, 8, 9]
display_model(model)
Image in a Jupyter notebook

StasFomin: Увы, моя вина, что нечетко (хотя однозначно, как в базе Вигго-Кана) сформулировал. Сейчас обновил формулировку, пояснив примером, о чем эта задача — т.е. это как покрытие, но с дополнительными ограничениями на различимость болезней. Увы, тут уже видно, что «9»-я болезнь не покрывается никакими тестами.

Модель в общем случае:


def MNT_ILP_model(initial_set, initial_subsets):
    """initial_set -- range(k) или range(a, b) со сдвигом"""
    model = p.ConcreteModel()
    
    model.InitialSet = initial_set
    model.SetLen = len(model.InitialSet)
    model.SetIndex = range(model.SetLen)
    
    model.SubsetCount = len(initial_subsets)
    model.SubsetIndex = range(model.SubsetCount)

    # model.InitialSubsets -- значения initial_subsets после сдвига 
    # range(a, b) -> range(b-a), для этого вычитаем из каждого элемента min.
    _min = min(model.InitialSet)
    model.InitialSubsets = [[elem - _min for elem in subset] for subset in initial_subsets]

    S = [[0 for _ in range(model.SetLen)] for _ in range(model.SubsetCount)]

    for i in range(model.SubsetCount):
        for elem in model.InitialSubsets[i]:
            S[i][elem] = 1

    model.Subsets = S

    model.x = p.Var(model.SubsetIndex, within=p.Binary)
    
    @model.Objective()
    def FinalCapacity(model):
        return sum(model.x[subset] for subset in model.SubsetIndex)

    @model.Constraint(model.SetIndex, model.SetIndex)
    def all_pair_of_elements_are_included_excluded(model, i, j):
        if i <= j:
            return p.Constraint.Feasible
        else:
            return sum(model.x[subset] * (model.Subsets[subset][i] ^ model.Subsets[subset][j]) for subset in model.SubsetIndex) >= 1

    return model

Трансляция к задаче о разрешимости

Если мощность найденного C′⊆C C' \subseteq C меньше, либо равна, заданного параметра γ(∣C∣)\gamma(|C|), то задача разрешима. Иначе нет.

def is_feasible_pyomo(initial_set, initial_subsets, gamma):
    model = MNT_ILP_model(initial_set, initial_subsets)

    solver = p.SolverFactory('cbc')
    try:
        solution = solver.solve(model)

        if str(solution['Solver'][0]['Termination condition']) == 'optimal' and model.FinalCapacity() <= gamma:
            return True, model
        else:
            return False, None
    except:
        return False, None

is_feasible_pyomo(gamma=100, initial_set=range(5), initial_subsets=[[0,1], [1]])[0]
WARNING: Loading a SolverResults object with a warning status into model.name="unknown"; - termination condition: infeasible - message from solver: <undefined>
False
is_feasible_pyomo(gamma = 4, initial_set=range(5), initial_subsets=[[0,1], [1], [2], [3], [4]])[0]
True
is_feasible_pyomo(gamma = 3, initial_set=range(5), initial_subsets=[[0,1], [1], [2], [3], [4]])[0]
False

Сведение 3SAT к задаче

Пусть нам дана формула Φ(x1,y1,z1,…,xp,yp,zp)=φ1(x1,y1,z1)⋀φ2(x2,y2,z2)…φp(xp,yp,zp)\Phi(x_1, y_1, z_1, \dots, x_p, y_p, z_p) = \varphi_1(x_1,y_1,z_1) \bigwedge \varphi_2(x_2,y_2,z_2) \dots \varphi_p(x_p,y_p,z_p) Здесь φ1(x1,y1,z1)=x1⋁y1⋁z1\varphi_1(x_1,y_1,z_1) = x_1 \bigvee y_1 \bigvee z_1 – конъюкция и x1,y1,z1x_1, y_1, z_1 – булевы литералы, либо их отрицание (например, φ1(x1,x4‾,x6‾)\varphi_1(x_1, \overline{x_4}, \overline{x_6})).

Пусть в Φ\Phi входит p конъюкций и k независимых булевых переменных (не считая их отрицаний).

Строим по 3SAT условие для MNT

Исходное множество значений C\textbf{C} для задачи Minimum Test Collection(MNT) будет состоять из следующих элементов:

  • "0" — вспомогательный нулевой элемент, не входит ни в одно подмножество семейства C\textbf{C},

  • "+∞+\infty" — вспомогательный элемент, содержится ровно в одном подмножестве C\textbf{C},

  • "{−1,−2,−3,…,−k}\{-1, -2, -3, \dots, -k\}" — "ядра" булевых переменных,

  • "{1,1.5,2,2.5,…,p,p+0.5}\{1, 1.5, 2, 2.5, \dots, p, p + 0.5 \}" — "ядра" конъюкций (целые) и вспомогательные элементы (дробные).

Строим исходное семейство подмножеств C\textbf{C}. В него будут входить следующие подмножества:

  • "черное" — содержит +∞+\infty и все положительные числа. N={+∞,1,1.5,…,p,p+0.5}∈C\textbf{N} = \{+\infty, 1, 1.5, \dots, p, p + 0.5\} \in \textbf{C},

  • "красные" — pp подмножеств, содержащих ровно 2 элемента: ii и i+0.5i + 0.5. R={r∣  r={i,i+0.5},i∈1…p‾}⊂C\textbf{R} = \Big\{ \textbf{r} |\; \textbf{r} = \big\{ i, i + 0.5 \big\}, i \in \overline{1\dots p} \Big\} \subset \textbf{C}

  • "синие" — 2k2k подмножеств, половина которых содержит число −i-i и conjiconj_i (назовем их xi\textbf{x}_\textbf{i}), а вторая половина — −i-i и conji‾\overline{conj_i} (назовем их xi‾\overline{\textbf{x}_\textbf{i}}). B={xi∣  xi={−i}∪conji,i∈1..k‾}∪{xi‾∣  xi‾={−i}∪conji‾,i∈1..k‾}⊂C\textbf{B} = \Big\{ \textbf{x}_\textbf{i} |\; \textbf{x}_\textbf{i} = \{ -i \} \cup conj_i, i \in \overline{1..k} \Big\} \cup \Big\{ \overline{\textbf{x}_\textbf{i}} |\; \overline{\textbf{x}_\textbf{i}} = \{ -i \} \cup \overline{conj_i}, i \in \overline{1..k} \Big\} \subset \textbf{C}

Утверждение:

Существование решения задачи Minimum Test Collection(MNT) с построенными S\textbf{S}, C\textbf{C} и gamma=1+p+kgamma = 1 + p + k эквивалентно существованию решения 3-SAT: Φ(x1,x1‾,x2,x2‾,…,xk,xk‾)=1\Phi(x_1, \overline{x_1}, x_2, \overline{x_2}, \dots, x_k, \overline{x_k}) = 1.

или
∃C′ – решение MNT(S,C,gamma)⟺∃{x1∗,x2∗,…xk∗} – решение 3SAT: Φ(x1,x1‾,x2,x2‾,…,xk,xk‾)=1\exists C' \text{ -- решение MNT}(S,C,gamma) \Longleftrightarrow \exists \{x^*_1, x^*_2, \dots x^*_k\} \text{ -- решение 3SAT: } \Phi(x_1, \overline{x_1}, x_2, \overline{x_2}, \dots, x_k, \overline{x_k}) = 1

Доказательство:

Влево: Пусть набор {x1∗,x2∗,…xk∗} – решение 3SAT: Φ(x1,x1‾,x2,x2‾,…,xk,xk‾)=1 \{x^*_1, x^*_2, \dots x^*_k\} \text{ -- решение 3SAT: } \Phi(x_1, \overline{x_1}, x_2, \overline{x_2}, \dots, x_k, \overline{x_k}) = 1 . Решим MNT с заданными S,C\textbf{S}, \textbf{C} и gamma=1+p+kgamma = 1 + p + k. Выберем из C\textbf{C} следующие подмножества:

  • "черное",

  • все "красные",

  • "синие" по следующему правилу: для каждого i∈1..k‾i \in \overline{1..k}, если xi∗=1x^*_i = 1 — берем {−i}∪conji∈B\{-i \} \cup conj_i \in \textbf{B}, если xi∗=0x^*_i = 0 — берем {−i}∪conji‾∈B\{-i \} \cup \overline{conj_i} \in \textbf{B}

Итого ∣C′∣=1+p+k |C'| = 1 + p + k подмножеств выбрано.

Докажем, что такой выбор подмножеств удовлетворяет условию MNT: Покажем, что каждый элемент из S\textbf{S} ни конфликтует ни с каким другим (т.е. для любой пары есть подмножество из |C'|, содержащее один элемент, но не второй):

  • Нулевой элемент

    • +∞+\infty: входит в "черное", а "0" — нет,

    • Отрицательные числа: ∀i∈1..k‾\forall i \in \overline{1..k} вы взяли "синее", содержащее −i-i, но не "0",

    • Положительные числа: каждое положительное число входит в "красное", куда не входит "0", "красные" взяли все.

  • +∞+\infty

    • Отрицательные числа: не входят в "черное", куда входит +∞+\infty,

    • Положительные числа: каждое положительное число входит в "красное", куда не входит +∞+\infty, "красные" взяли все.

  • Отрицательные числа

    • Другие отрицательные: ∀i∈1..k‾\forall i \in \overline{1..k} вы взяли "синее", содержащее единственное отрицательное число: −i-i,

    • Положительные числа: входят в "черное", куда не входят отрицательные.

  • Дробные положительные:

    • Другие дробные: в каждом "красном" ровно одно дробное число, "красные" взяли все.

    • Целые положительные (не из "красного", в которое входит исходное дробное): берем это "красное".

  • Целые положительные:

    • Другие целые положительные: в каждом "красном" ровно одно целое число, "красные" взяли все.

Итак, не разобрана только одна связь: дробные положительные и целые положительные, входящие в одно "красное" подмножество. Пусть в 3SAT-задаче φi=xa⋁xb⋁xc\varphi_i = x_a \bigvee x_b \bigvee x_c. Тогда, по построению, целое число ii, помимо 1 "черного" и 1 "красного" подмножества, потенциально входит в 3 "синих": в те, что содержат числа −a,−b,−c-a, -b, -c ({−a}∪conja\{-a\} \cup conj_a и т.д.). При этом, исходная задача разрешима, значит xa∗=1x^*_a = 1, или xb∗=1x^*_b = 1, или xc∗=1x^*_c = 1, значит, по построению "синих", мы выбрали в C′C' одно из этих трёх "синих". Таким образом это "синее" подмножество разрешает конфликт между ii и i+0.5i + 0.5.

Замечание: если φi=xa‾⋁xb⋁xc‾\varphi_i = \overline{x_a} \bigvee x_b \bigvee \overline{x_c}, идея выше не поменяется, в таком случае xa∗=0x^*_a = 0, или xb∗=1x^*_b = 1, или xc∗=0x^*_c = 0 и ii входит в одно из ({−a}∪conja‾\{-a\} \cup \overline{conj_a}, или {−b}∪conjb\{-b\} \cup conj_b, или {−c}∪conjc‾\{-c\} \cup \overline{conj_c}).

Вправо:

Пусть C′C' — решение MNT( S\mathbf{S}, C\mathbf{C}, gamma=1+p+kgamma = 1 + p + k). Заметим, что исходя из построенного S\mathbf{S} и C\mathbf{C}, чтобы не было конфликтов между элементами S\mathbf{S}:

  1. "Черное" подмножество входит в C′C', иначе конфликт между "0" и +∞+\infty ("0" не входит ни в одно, +∞+\infty входит только в "черное" ).

  2. Чтобы исключить конфликт между +∞+\infty и дробными положительными, а так как каждое дробного положительного входит только в "черное" (что не решает конфликт) и в одно "красное", необходимо внести все "красные" в C′C'. Мы уже потратили 1+p1 + p подмножеств.

  3. Чтобы исключить конфликт между "0" и отрицательными числами, для каждого i∈1..k‾i \in \overline{1..k} необходимо взять в C′C' либо xi={−i}∪conji\textbf{x}_\textbf{i} = \{-i\} \cup conj_i, либо xi‾={−i}∪conji‾\overline{\textbf{x}_\textbf{i}} = \{-i\} \cup \overline{conj_i}. При этом нам осталось доступно только gamma−1−p=kgamma - 1 - p = k подмножеств, значит для каждого i∈1..k‾i \in \overline{1..k} надо взять одно и ровно одно: xi\textbf{x}_\textbf{i} или xi‾\overline{\textbf{x}_\textbf{i}}.

  4. Чтобы устранить конфликт между дробными и целыми положительными, находящимися в одной "красной" (т.е. между ii и i+0.5i + 0.5), для каждого i∈1..p‾i \in \overline{1..p} в C′C' должно быть хотя бы одно "синее" подмножество, содержащее ii.

Докажем, что из этих выводов следует, что набор {x1∗,x2∗,…xk∗}\{x^*_1, x^*_2, \dots x^*_k\}, где xi∗={1,мы выбрали xi(xi∈C′)0,мы выбрали xi‾(xi‾∈C′)\begin{equation} x^*_i = \begin{cases} 1 &, \text{мы выбрали }\textbf{x}_\textbf{i} (\textbf{x}_\textbf{i} \in C' )\\ 0 &, \text{мы выбрали } \overline{\textbf{x}_\textbf{i}} (\overline{\textbf{x}_\textbf{i}} \in C') \end{cases} \end{equation} является решением исходной 3SAT.

Из (3) следует, что xi∗x^*_i определяется однозначно ∀i∈1..k‾\forall i \in \overline{1..k}. Из (4) следует, что ∀i∈1..p‾  ∃\forall i \in \overline{1..p}\; \exists хотя бы одно "синее" подмножество, содержащее ii, лежащее в C′C'. Пусть это синее подмножество сформировано вокруг элемента −k_i-k\_i, то есть это xk_i\textbf{x}_\textbf{k\_i}, если xk_ix_{k\_i} входит в φi\varphi_i, либо xk_i‾\overline{\textbf{x}_\textbf{k\_i}}, если в φi\varphi_i входит xk_i‾\overline{x_{k\_i}} . Таким образом: [∀i∈1..p‾  ∃xk_i∈C′:  xk_i∈φi∀i∈1..p‾  ∃xk_i‾∈C′:  xk_i‾∈φi \left[ \begin{gathered} \forall i \in \overline{1..p} \; \exists \textbf{x}_\textbf{k\_i} \in C':\; x_{k\_i} \in \varphi_i \\ \forall i \in \overline{1..p} \; \exists \overline{\textbf{x}_\textbf{k\_i}} \in C':\; \overline{x_{k\_i}} \in \varphi_i \end{gathered} \right.

То есть для каждой конъюкции φi\varphi_i из Φ=φ1∧φ2∧⋯∧φp\Phi = \varphi_1 \wedge \varphi_2 \wedge \dots \wedge \varphi_p (т.е. для каждого положительного целого элемента) есть хотя бы один элемент xk_i∗=1x^*_{k\_i} = 1 (либо xk_i∗‾=0‾=1\overline{x^*_{k\_i}} = \overline{0} = 1), входящий в эту конъюкцию, значит все конъюкции выполняются, значит {x1∗,x2∗,…xk∗}\{x^*_1, x^*_2, \dots x^*_k\} — решение Φ=1\Phi = 1.

Визуализируем

Функция, генерирующая рандомную формулу 3-КНФ 

def Generate_random_3cnf(seed, conjunction_count, variables_count):
    random.seed(seed)
    return CNF(from_clauses=psc.rand3cnf(conjunction_count, variables_count))

cnf = Generate_random_3cnf(3, 4, 5)
cnf.clauses
[(-2, 5, 4), (-3, 1, 5), (3, 5, -4), (5, 3, -1)]

Функция, превращающая любую 3-КНФ формулу в экземляр задачи MNT (тройку (S, C, gamma)) по алгоритму, описанному выше.

def transform_3cnf_to_MNT(cnf3):
    clauses = cnf3.clauses

    var_names = set()
    for clause in clauses:
        for var in clause:
            var_names.add(abs(var))
    # psc.rand3cnf всегда называет переменные в порядке от 1 до k (без пробелов)
    
    # Число переменных k
    k = len(var_names)
    # Число конъюкций p
    p = len(clauses)
    # Параметр gamma
    gamma = 1 + k + p

    # Исходное множество значений S.
    # Для простоты положительные числа умножены на 2 и уменьшины на 1 (теперь они не range(1, p + 0.5, step=0.5), а range(1, 2p + 1)),
    # роль infty выполняет 2p + 1.
    S = list(range(-k, 0)) + [0] + list(range(1, 2 * p + 1 )) + [2 * p + 1]
    # Исходное семейство подмножеств С.
    C = []
    # Добавляем "черное"
    C.append(list(range(1, 2 * p + 2)))
    # Добавляем "красные"
    for i in range(0, p):
        C.append([2 * i + 1, 2 * (i + 1)])
    # Добавляем "синие"
    var_in_conj, non_var_in_conj = [[-i] for i in range(1, k + 1)], [[-i] for i in range(1, k + 1)]
    for conj_id, clause in enumerate(clauses):
        for var in clause:
            if var > 0:
                # var_in_conj[0] = [-1]
                var_in_conj[abs(var) - 1].append(2 * conj_id + 1)
            else:
                non_var_in_conj[abs(var) - 1].append(2 * conj_id + 1)
    C.extend(var_in_conj)
    C.extend(non_var_in_conj)
    return S, C, gamma

Пример того, что выдает эта функция.

cnf3 = Generate_random_3cnf(seed=1, conjunction_count=3, variables_count=2)

print(f'Исходная задача: {cnf3.clauses}')

(initial_set, initial_subsets, gamma) = transform_3cnf_to_MNT(cnf3)

print(f'initial_set: {initial_set}')
print(f'initial_subsets: {initial_subsets}')
print(f'gamma: {gamma}')

(is_feasible, model) = is_feasible_pyomo(
    initial_set, initial_subsets, gamma
)
Исходная задача: [(-1, -2), (1, -2), (2, -1)] initial_set: [-2, -1, 0, 1, 2, 3, 4, 5, 6, 7] initial_subsets: [[1, 2, 3, 4, 5, 6, 7], [1, 2], [3, 4], [5, 6], [-1, 3], [-2, 5], [-1, 1, 5], [-2, 1, 3]] gamma: 6 
if is_feasible:
    selected_subsets = [initial_subsets[subset] for subset in model.SubsetIndex if round(model.x[subset]()) == 1]
    print(f'selected_subsets: {selected_subsets}')
    svg = display_model(model)
svg
selected_subsets: [[1, 2, 3, 4, 5, 6, 7], [1, 2], [3, 4], [5, 6], [-1, 1, 5], [-2, 1, 3]]
Image in a Jupyter notebook

Проверка выполнимости через ЦЛП и SAT-Solver

Проверка разрешимости 3-КНФ через преобразование к MNT и решение через pyomo-ЦЛП

def is_feasible_ILP(cnf3):
    (initial_set, initial_subsets, gamma) = transform_3cnf_to_MNT(cnf3)
    return is_feasible_pyomo(initial_set, initial_subsets, gamma)[0]

Проверка разрешимости 3-КНФ через pysat Solver

def is_feasible_SAT(cnf3):
    from pysat.solvers import Solver

    solver = Solver(bootstrap_with=cnf3)
    return solver.solve()

Проверим что это вообще работает на тривиальных 3-КНФ.

true_cnf = CNF(from_clauses=[(1,2,3), (1,-2,-3), (1,2,-3)])
false_cnf = CNF(from_clauses=[(1,2,3), (1,2,-3), (1,-2, 3), (1,-2,-3), (-1,2,3), (-1,2,-3), (-1,-2, 3), (-1,-2,-3),])

is_feasible_SAT(true_cnf), is_feasible_SAT(false_cnf), is_feasible_ILP(true_cnf), is_feasible_ILP(false_cnf),
(True, False, True, False)

Тест, который генерирует test_lengthtest\_length случайных 3-КНФ формул длины до max_conj_countmax\_conj\_count и проверяет, что наш ILP и pysat солвер выдают один и тот же ответ на задачу разрешимости


def run_test(test_length, max_conj_count):
    start = time.time()
    misses_count = 0
    for _ in range(test_length):
        seed = randint(0, 10_000_000)
        conj_count = randint(1, max_conj_count)
        var_count = randint(1, max(1, conj_count // 3))
        Phi = Generate_random_3cnf(seed=seed, conjunction_count=conj_count, variables_count=var_count)

        if is_feasible_ILP(Phi) != is_feasible_SAT(Phi):
            misses_count += 1
            if misses_count % 5 == 0:
                print("Got 5 misses")
    return misses_count, time.time() - start

Теперь просто гоняем тесты

Тест на больших формул:

run_test(10, 50)
(0, 12.412810325622559)
run_test(120, 50)
(0, 156.25865602493286)

Много тестов для формулы поменьше:

run_test(500, 10)
(0, 100.30018091201782)

Судя по результатам, для всех прогнанных тестов задачи были эквивалентны.

Теперь хочется прогнать на совсем больших формулах (и посмотреть, что задачка действительно NPNP: с увеличением входа время растет экспоненциально).

def big_test(conj_count):
    mark1 = time.time()

    Phi = Generate_random_3cnf(seed=randint(0,1_000_000), conjunction_count=conj_count, variables_count=conj_count)
    mark2 = time.time()

    ILP_res = is_feasible_ILP(Phi)
    mark3 = time.time()

    SAT_res = is_feasible_SAT(Phi)
    mark4 = time.time()

    if ILP_res != SAT_res:
        print("Err")
    else:
        print("Ok")
    return mark2 - mark1, mark3 - mark2, mark4 - mark3

big_test(50)
Ok
(0.005578756332397461, 11.565046072006226, 0.0005283355712890625)
big_test(100)
Ok
(0.06852245330810547, 77.0499587059021, 0.00041866302490234375)
big_test(150)
Ok
(0.012944936752319336, 262.47025990486145, 0.0013048648834228516)
# big_test(200)

Судя по времени выполнения, решение 3SAT задачки через ЦЛП — не лучшая затея.